/********************************************************
 * @file: vector.h
 * @date: 25/09/2008
 ********************************************************
 * FBElem v1_0 - (c)copyright 2008 GoBEM
 ********************************************************/

#ifndef _VECTOR_H_
#define _VECTOR_H_

// includes
//----------
#include <iostream>
#include <assert.h>
#include "matrix.h"
#include <math.h>

// forward declarations
//----------------------
template <class T> 
class Matrix;

// useful typedefs
//-----------------
using namespace std;


/*!
 * CLASS
 * Vector
 * 
 * This class defines the capabilities of a vector (a one dimensional array).
 *
 * It has several functions which aid the vector manipulations and operations,
 * is an efficient way. But on should think in use this as a wrapper class in
 * order to have access to different LIBs.
 *
 * Since it is an template, it allows the use of both double and complex matrices.
 * 
 * @author Andre Pereira (andremaues@gmail.com)
 * @version 1.0.0
 * @see Vector
 */
template <class T>
class Vector
{

public:

  // Constructors
  Vector( unsigned int _length = 2);
  Vector(const Vector<T>& _v);

  // Destructor
  ~Vector();

  // Assignment operator
	Vector<T>& operator = (const T& _a);
	Vector<T>& operator = (const Vector<T>& _v);

  // Value extraction method
  unsigned int getLength() const;

  // Subscript operator
  T& operator () (unsigned int _index);
  const T& operator () (unsigned int _index) const;

  // Calculation operators
  Vector<T> operator + (const Vector<T>& _v) const;
  Vector<T> operator - (const Vector<T>& _v) const;
  Matrix<T> operator * (const Matrix<T>& _m) const;
  Vector<T> operator * (const T& _a) const;
  Vector<T> operator / (const T& _a) const;
	T operator * (const Vector<T>& _v) const;

  // Combined assignment - calculation operators
  Vector<T>& operator += (const Vector<T>& _v);
  Vector<T>& operator -= (const Vector<T>& _v);
  Vector<T>& operator *= (const T& _a);

  // Equality operators
	bool operator == (const Vector<T>& _v) const;

  // Miscellaneous -methods
  void zero();
  //void unit();

  // Utility methods
	void resize(unsigned int _n);
	double norm() const;
	Vector<T> cross( const Vector<T>& _v ) const;
         T maxV() const;
        T minV() const;
	void print() const;
  // norm(), dot()

  // Type of vectors


private:

  T * m_data;
	unsigned int m_length;

};

typedef Vector<double> DVector;
typedef Vector<int>    IVector;


//-----------------------------------------------------
template <class T> inline
Vector<T>::Vector( unsigned int _length )
{
	m_data = new T[ _length ];
	assert( m_data != 0 );
	m_length = _length;
}
//-----------------------------------------------------

//-----------------------------------------------------
template <class T> inline
Vector<T>::Vector(const Vector<T>& v)
{
	m_data = new T[m_length=v.m_length];
	assert(m_data!=0);
	for(unsigned int i=0;i<m_length;++i)
    m_data[i]=v.m_data[i];
}
//-----------------------------------------------------

//-----------------------------------------------------
template <class T> inline
Vector<T>::~Vector()
{
  delete [] m_data;
  m_data = 0;
}
//-----------------------------------------------------

//-----------------------------------------------------
template <class T> inline
Vector<T>&
Vector<T>::operator = (const T& _a)
{
  for( unsigned int i=0; i<m_length; ++i )
    m_data[i] = _a;
	return *this;
}
//-----------------------------------------------------

//-----------------------------------------------------
template <class T> inline
Vector<T>&
Vector<T>::operator = (const Vector<T>& v)
{
	if(this != &v)
	{
		delete [] m_data;
		m_data= new T[m_length=v.m_length];
		assert(m_data!=0);
		for(unsigned int i=0;i<m_length;++i)
      m_data[i]=v.m_data[i];
	}
	return *this;
}
//-----------------------------------------------------

//-----------------------------------------------------
template <class T> inline
unsigned int
Vector<T>::getLength() const
{
  return m_length;
}
//-----------------------------------------------------

//-----------------------------------------------------
template <class T> inline
T&
Vector<T>::operator () (unsigned int _index )
{
	assert( _index < m_length );
  return m_data[_index];
}
//-----------------------------------------------------

//-----------------------------------------------------
template <class T> inline
const T&
Vector<T>::operator () (unsigned int _index ) const
{
	assert( _index < m_length );
  return m_data[_index];
}
//-----------------------------------------------------

//-----------------------------------------------------
template <class T> inline
Vector<T>
Vector<T>::operator + ( const Vector<T>& v ) const
{
	assert( v.m_length == m_length );
  Vector<T> sum( m_length );
	for(unsigned int i=0; i<m_length ;++i)
    sum.m_data[i] = m_data[i] + v.m_data[i];
  return sum;
}
//-----------------------------------------------------

//-----------------------------------------------------
template <class T> inline
Vector<T>
Vector<T>::operator - ( const Vector<T>& v ) const
{
	assert( v.m_length == m_length );
  Vector<T> dif( m_length );
	for(unsigned int i=0; i<m_length ;++i)
    dif.m_data[i] = m_data[i] - v.m_data[i];
  return dif;
}
//-----------------------------------------------------

//-----------------------------------------------------
template <class T> inline
Matrix<T>
Vector<T>::operator * ( const Matrix<T>& m ) const
{
	assert( m.getRows() == 1 );
  unsigned int C = m.getCols();
  Matrix<T> mult( m_length, C );
  unsigned int i, j;
  for( i=0; i<m_length; i++ )
	  for( j=0; j<C; j++ )
      mult(i,j) = m_data[i] * m(0,j);
  return mult;
}
//-----------------------------------------------------

//-----------------------------------------------------
template <class T> inline
Vector<T>
Vector<T>::operator * ( const T& a ) const
{
  Vector<T> result(m_length);
  for( unsigned int i=0;i<m_length;i++ )
    result.m_data[i] = m_data[i]*a;
  return result;
}
//-----------------------------------------------------

//-----------------------------------------------------
template <class T> inline
Vector<T>
Vector<T>::operator / ( const T& a ) const
{
  Vector<T> result(m_length);
  for( unsigned int i=0;i<m_length;i++ )
    result.m_data[i] = m_data[i]/a;
  return result;
}
//-----------------------------------------------------

//-----------------------------------------------------
template <class T> inline
T
Vector<T>::operator * ( const Vector<T>& v ) const
{
	assert( v.m_length == m_length );
  T dot = 0.0;
	for(unsigned int i=0; i<m_length ;++i)
    dot += m_data[i] * v.m_data[i];
  return dot;
}
//-----------------------------------------------------

//-----------------------------------------------------
template <class T>
Vector<T>&
Vector<T>::operator += (const Vector<T> &v)
{
	assert( v.m_length == m_length );
  for( unsigned int i=0;i<m_length; i++)
    m_data[i] += v.m_data[i];
  return *this;
}
//-----------------------------------------------------
     
//-----------------------------------------------------
template <class T>
Vector<T>& 
Vector<T>::operator -= (const Vector<T> &v)
{
	assert( v.m_length == m_length );
  for( unsigned int i=0;i<m_length; i++)
    m_data[i] -= v.m_data[i];
  return *this;
}
//-----------------------------------------------------

//-----------------------------------------------------
template <class T>
Vector<T>& 
Vector<T>::operator *= (const T& a)
{
  for( unsigned int i=0;i<m_length;i++ )
    m_data[i] *= a;
  return *this;
}
//-----------------------------------------------------

//-----------------------------------------------------
template <class T>
bool
Vector<T>::operator == (const Vector<T>& v) const
{
	if(m_length != v.m_length) return false;
	for(unsigned int i=0;i<m_length;i++)
    if(m_data[i] != v.m_data[i])
      return false;
	return true;
}
//-----------------------------------------------------

//-----------------------------------------------------
template <class T> inline
void
Vector<T>::zero()
{
  for( unsigned int i=0; i<m_length; ++i )
    m_data[i] = T(0);
}
//-----------------------------------------------------

//-----------------------------------------------------
template <class T> inline
void
Vector<T>::resize(unsigned int _n)
{
	delete [] m_data;
	m_data = new T[m_length=_n];
	assert(m_data !=0);
}
//-----------------------------------------------------

//-----------------------------------------------------
template <class T> inline
double
Vector<T>::norm() const
{
  return sqrt( *this * *this );
}
//-----------------------------------------------------

//-----------------------------------------------------
template <class T> inline
Vector<T>
Vector<T>::cross( const Vector<T>& _v ) const
{
  assert( _v.m_length == m_length );
  assert( m_length <= 3 && m_length > 1 );
  Vector<T> crossProd(3);
  if( m_length == 3 )
  {
    crossProd.m_data[0] = m_data[1]*_v.m_data[2] - m_data[2]*_v.m_data[1];
    crossProd.m_data[1] = - m_data[0]*_v.m_data[2] - m_data[2]*_v.m_data[0];
  }
  else
  {
    crossProd.m_data[0] = 0.0;
    crossProd.m_data[1] = 0.0;
  }
  crossProd.m_data[2] = m_data[0]*_v.m_data[1] - m_data[1]*_v.m_data[0];
  return crossProd;
}
//-----------------------------------------------------

//-----------------------------------------------------
template <class T> inline
T
Vector<T>::maxV() const
{
  T maxVal = m_data[0];
  for(unsigned int i=0;i<m_length;i++)
    if(maxVal < m_data[i])
      maxVal = m_data[i];
  return maxVal;
}
//-----------------------------------------------------

//-----------------------------------------------------
template <class T> inline
T
Vector<T>::minV() const
{
  T minVal = m_data[0];
  for(unsigned int i=0;i<m_length;i++)
    if(minVal > m_data[i])
      minVal = m_data[i];
  return minVal;
}
//-----------------------------------------------------

//-----------------------------------------------------
template <class T> inline
void
Vector<T>::print() const
{
  unsigned int i;
  T * p;
	for(i=0, p=m_data; i<m_length; ++i)
    cout << *p++ << '\n';
	cout << endl;
}
//-----------------------------------------------------

////////////////////////////////////////////////
// GLOBAL OPERATORS
////////////////////////////////////////////////

//-----------------------------------------------------
template <class T> inline
ostream&
operator << (ostream& s, const Vector<T>& v)
{
  for( unsigned int i=0; i < v.getLength(); ++i )
    s << v(i) << '\t';
  s << endl;
  return s;
}
//-----------------------------------------------------

//-----------------------------------------------------
template <class T> inline
Vector<T>
operator * (const T &a, const Vector<T> &x)
{
  int N = x.getLength();
  Vector<T> result(N);
  for (int i=0;i<N;i++)
    result(i) = x(i)*a;
  return result;
}
//-----------------------------------------------------
/*
//-----------------------------------------------------
template <class T> inline
Vector<T>&
operator = (Vector<T>& u, const Vector<T>& v)
{
	if(&u != &v)
	{
    unsigned int L = v.getLength();
    u.resize( L )
		for(unsigned int i=0;i<L;++i)
      u(i) = v(i);
	}
	return u;
}
//-----------------------------------------------------
*/
/*
          
template <class TYPE>
TYPE dot(const MV_Vector<TYPE> &x, const MV_Vector<TYPE> &y)
{
        
  //  Check for compatible dimensions:
  if (x.size() != y.size())
      {
         cout << "Incompatible dimensions in dot(). " << endl;
         exit(1);
      }

      TYPE temp=0.0;
      for (int i=0; i<x.size();i++)
           temp += x(i)*y(i);
      return temp;
}

template <class TYPE>
TYPE norm(const MV_Vector<TYPE> &x)
{
      TYPE temp = dot(x,x);
      return sqrt(temp);
}
*/

#endif // _MATRIX_H_
